Setting Boundaries for Statistical Mechanics.

Statistical mechanics has grown without bounds in space. Statistical mechanics of noninteracting point particles in an unbounded perfect gas is widely used to describe liquids like concentrated salt solutions of life and electrochemical technology, including batteries. Liquids are filled with interacting molecules. A perfect gas is a poor model of a liquid. Statistical mechanics without spatial bounds is impossible as well as imperfect, if molecules interact as charged particles, as nearly all atoms do. The behavior of charged particles is not defined until boundary structures and values are defined because charges are governed by Maxwell's partial differential equations. Partial differential equations require boundary structures and conditions. Boundary conditions cannot be defined uniquely 'at infinity' because the limiting process that defines 'infinity' includes such a wide variety of structures and behaviors, from elongated ellipses to circles, from light waves that never decay, to dipolar fields that decay steeply, to Coulomb fields that hardly decay at all. Boundaries and boundary conditions needed to describe matter are not prominent in classical statistical mechanics. Statistical mechanics of bounded systems is described in the EnVarA system of variational mechanics developed by Chun Liu, more than anyone else. EnVarA treatment does not yet include Maxwell equations.

Effects of Diffusion Coefficients and Permanent Charge on Reversal Potentials in Ionic Channels.

In this work, the dependence of reversal potentials and zero-current fluxes on diffusion coefficients are examined for ionic flows through membrane channels. The study is conducted for the setup of a simple structure defined by the profile of permanent charges with two mobile ion species, one positively charged (cation) and one negatively charged (anion). Numerical observations are obtained from analytical results established using geometric singular perturbation analysis of classical Poisson-Nernst-Planck models. For 1:1 ionic mixtures with arbitrary diffusion constants, Mofidi and Liu (arXiv:1909.01192) conducted a rigorous mathematical analysis and derived an equation for reversal potentials. We summarize and extend these results with numerical observations for biological relevant situations. The numerical investigations on profiles of the electrochemical potentials, ion concentrations, and electrical potential across ion channels are also presented for the zero-current case. Moreover, the dependence of current and fluxes on voltages and permanent charges is investigated. In the opinion of the authors, many results in the paper are not intuitive, and it is difficult, if not impossible, to reveal all cases without investigations of this type.

Molecular Mean-Field Theory of Ionic Solutions: A Poisson-Nernst-Planck-Bikerman Model.

We have developed a molecular mean-field theory-fourth-order Poisson-Nernst-Planck-Bikerman theory-for modeling ionic and water flows in biological ion channels by treating ions and water molecules of any volume and shape with interstitial voids, polarization of water, and ion-ion and ion-water correlations. The theory can also be used to study thermodynamic and electrokinetic properties of electrolyte solutions in batteries, fuel cells, nanopores, porous media including cement, geothermal brines, the oceanic system, etc. The theory can compute electric and steric energies from all atoms in a protein and all ions and water molecules in a channel pore while keeping electrolyte solutions in the extra- and intracellular baths as a continuum dielectric medium with complex properties that mimic experimental data. The theory has been verified with experiments and molecular dynamics data from the gramicidin A channel, L-type calcium channel, potassium channel, and sodium/calcium exchanger with real structures from the Protein Data Bank. It was also verified with the experimental or Monte Carlo data of electric double-layer differential capacitance and ion activities in aqueous electrolyte solutions. We give an in-depth review of the literature about the most novel properties of the theory, namely Fermi distributions of water and ions as classical particles with excluded volumes and dynamic correlations that depend on salt concentration, composition, temperature, pressure, far-field boundary conditions etc. in a complex and complicated way as reported in a wide range of experiments. The dynamic correlations are self-consistent output functions from a fourth-order differential operator that describes ion-ion and ion-water correlations, the dielectric response (permittivity) of ionic solutions, and the polarization of water molecules with a single correlation length parameter.

Field theory of reaction-diffusion: Law of mass action with an energetic variational approach.

We extend the energetic variational approach so it can be applied to a chemical reaction system with general mass action kinetics. Our approach starts with an energy-dissipation law. We show that the chemical equilibrium is determined by the choice of the free energy and the dynamics of the chemical reaction is determined by the choice of the dissipation. This approach enables us to couple chemical reactions with other effects, such as diffusion and drift in an electric field. As an illustration, we apply our approach to a nonequilibrium reaction-diffusion system in a specific but canonical setup. We show by numerical simulations that the input-output relation of such a system depends on the choice of the dissipation.

Do Bistable Steric Poisson-Nernst-Planck Models Describe Single-Channel Gating?

Experiments measuring currents through single protein channels show unstable currents, a phenomena called the gating of a single channel. Channels switch between an "open" state with a well-defined single amplitude of current and "closed" states with nearly zero current. The existing mean-field theory of ion channels focuses almost solely on the open state. The physical modeling of the dynamical features of ion channels is still in its infancy and does not describe the transitions between open and closed states nor the distribution of the duration times of open states. One hypothesis is that gating corresponds to noise-induced fast transitions between multiple steady (equilibrium) states of the underlying system. In this work, we aim to test this hypothesis. Particularly, our study focuses on the (high-order) steric Poisson-Nernst-Planck (PNP)-Cahn-Hilliard model since it has been successful in predicting permeability and selectivity of ionic channels in their open state and since it gives rise to multiple steady states. We show that this system gives rise to a gatinglike behavior, but that important features of this switching behavior are different from the defining features of gating in biological systems. Furthermore, we show that noise prohibits switching in the system of study. The above phenomena are expected to occur in other PNP-type models, strongly suggesting that one has to go beyond overdamped (gradient flow) Nernst-Planck type dynamics to explain the spontaneous gating of single channels.

Poisson-Fermi modeling of ion activities in aqueous single and mixed electrolyte solutions at variable temperature.

The combinatorial explosion of empirical parameters in tens of thousands presents a tremendous challenge for extended Debye-Hückel models to calculate activity coefficients of aqueous mixtures of the most important salts in chemistry. The explosion of parameters originates from the phenomenological extension of the Debye-Hückel theory that does not take steric and correlation effects of ions and water into account. By contrast, the Poisson-Fermi theory developed in recent years treats ions and water molecules as nonuniform hard spheres of any size with interstitial voids and includes ion-water and ion-ion correlations. We present a Poisson-Fermi model and numerical methods for calculating the individual or mean activity coefficient of electrolyte solutions with any arbitrary number of ionic species in a large range of salt concentrations and temperatures. For each activity-concentration curve, we show that the Poisson-Fermi model requires only three unchanging parameters at most to well fit the corresponding experimental data. The three parameters are associated with the Born radius of the solvation energy of an ion in electrolyte solution that changes with salt concentrations in a highly nonlinear manner.

Nonlocal Poisson-Fermi model for ionic solvent.

We propose a nonlocal Poisson-Fermi model for ionic solvent that includes ion size effects and polarization correlations among water molecules in the calculation of electrostatic potential. It includes the previous Poisson-Fermi models as special cases, and its solution is the convolution of a solution of the corresponding nonlocal Poisson dielectric model with a Yukawa-like kernel function. The Fermi distribution is shown to be a set of optimal ionic concentration functions in the sense of minimizing an electrostatic potential free energy. Numerical results are reported to show the difference between a Poisson-Fermi solution and a corresponding Poisson solution.

Poisson-Fermi Modeling of the Ion Exchange Mechanism of the Sodium/Calcium Exchanger.

The ion exchange mechanism of the sodium/calcium exchanger (NCX) crystallized by Liao et al. in 2012 is studied using the Poisson-Fermi theory developed by Liu and Eisenberg in 2014. A cycle of binding and unbinding is proposed to account for the Na(+)/Ca(2+) exchange function of the NCX molecule. Outputs of the theory include electric and steric fields of ions with different sizes, correlations of ions of different charges, and polarization of water, along with number densities of ions, water molecules, and interstitial voids. We calculate the electrostatic and steric potentials of the four binding sites in NCX, i.e., three Na(+) binding sites and one Ca(2+) binding site, with protein charges provided by the software PDB2PQR. The energy profiles of Na(+) and Ca(2+) ions along their respective Na(+) and Ca(2+) pathways in experimental conditions enable us to explain the fundamental mechanism of NCX that extrudes intracellular Ca(2+) across the cell membrane against its chemical gradient by using the downhill gradient of Na(+). Atomic and numerical details of the binding sites are given to illustrate the 3 Na(+):1 Ca(2+) stoichiometry of NCX. The protein NCX is a catalyst. It does not provide (free) energy for transport. All energy for transport in our model comes from the ions in surrounding baths.

Numerical methods for a Poisson-Nernst-Planck-Fermi model of biological ion channels.

Numerical methods are proposed for an advanced Poisson-Nernst-Planck-Fermi (PNPF) model for studying ion transport through biological ion channels. PNPF contains many more correlations than most models and simulations of channels, because it includes water and calculates dielectric properties consistently as outputs. This model accounts for the steric effect of ions and water molecules with different sizes and interstitial voids, the correlation effect of crowded ions with different valences, and the screening effect of polarized water molecules in an inhomogeneous aqueous electrolyte. The steric energy is shown to be comparable to the electrical energy under physiological conditions, demonstrating the crucial role of the excluded volume of particles and the voids in the natural function of channel proteins. Water is shown to play a critical role in both correlation and steric effects in the model. We extend the classical Scharfetter-Gummel (SG) method for semiconductor devices to include the steric potential for ion channels, which is a fundamental physical property not present in semiconductors. Together with a simplified matched interface and boundary (SMIB) method for treating molecular surfaces and singular charges of channel proteins, the extended SG method is shown to exhibit important features in flow simulations such as optimal convergence, efficient nonlinear iterations, and physical conservation. The generalized SG stability condition shows why the standard discretization (without SG exponential fitting) of NP equations may fail and that divalent Ca(2+) may cause more unstable discrete Ca(2+) fluxes than that of monovalent Na(+). Two different methods-called the SMIB and multiscale methods-are proposed for two different types of channels, namely, the gramicidin A channel and an L-type calcium channel, depending on whether water is allowed to pass through the channel. Numerical methods are first validated with constructed models whose exact solutions are known. The experimental data of both channels are then used to verify and explain novel features of PNPF as compared with previous PNP models. The PNPF currents are in accord with the experimental I-V (V for applied voltages) data of the gramicidin A channel and I-C (C for bath concentrations) data of the calcium channel with 10(-8)-fold bath concentrations that pose severe challenges in theoretical simulations.

Poisson-Nernst-Planck-Fermi theory for modeling biological ion channels.

A Poisson-Nernst-Planck-Fermi (PNPF) theory is developed for studying ionic transport through biological ion channels. Our goal is to deal with the finite size of particle using a Fermi like distribution without calculating the forces between the particles, because they are both expensive and tricky to compute. We include the steric effect of ions and water molecules with nonuniform sizes and interstitial voids, the correlation effect of crowded ions with different valences, and the screening effect of water molecules in an inhomogeneous aqueous electrolyte. Including the finite volume of water and the voids between particles is an important new part of the theory presented here. Fermi like distributions of all particle species are derived from the volume exclusion of classical particles. Volume exclusion and the resulting saturation phenomena are especially important to describe the binding and permeation mechanisms of ions in a narrow channel pore. The Gibbs free energy of the Fermi distribution reduces to that of a Boltzmann distribution when these effects are not considered. The classical Gibbs entropy is extended to a new entropy form - called Gibbs-Fermi entropy - that describes mixing configurations of all finite size particles and voids in a thermodynamic system where microstates do not have equal probabilities. The PNPF model describes the dynamic flow of ions, water molecules, as well as voids with electric fields and protein charges. The model also provides a quantitative mean-field description of the charge/space competition mechanism of particles within the highly charged and crowded channel pore. The PNPF results are in good accord with experimental currents recorded in a 10(8)-fold range of Ca(2+) concentrations. The results illustrate the anomalous mole fraction effect, a signature of L-type calcium channels. Moreover, numerical results concerning water density, dielectric permittivity, void volume, and steric energy provide useful details to study a variety of physical mechanisms ranging from binding, to permeation, blocking, flexibility, and charge/space competition of the channel.

Analytical models of calcium binding in a calcium channel.

The anomalous mole fraction effect of L-type calcium channels is analyzed using a Fermi like distribution with the experimental data of Almers and McCleskey [J. Physiol. 353, 585 (1984)] and the atomic resolution model of Lipkind and Fozzard [Biochemistry 40, 6786 (2001)] of the selectivity filter of the channel. Much of the analysis is algebraic, independent of differential equations. The Fermi distribution is derived from the configuration entropy of ions and water molecules with different sizes, different valences, and interstitial voids between particles. It allows us to calculate potentials and distances (between the binding ion and the oxygen ions of the glutamate side chains) directly from the experimental data using algebraic formulas. The spatial resolution of these results is comparable with those of molecular models, but of course the accuracy is no better than that implied by the experimental data. The glutamate side chains in our model are flexible enough to accommodate different types of binding ions in different bath conditions. The binding curves of Na(+) and Ca(2+) for [CaCl2] ranging from 10(-8) to 10(-2) M with a fixed 32 mM background [NaCl] are shown to agree with published Monte Carlo simulations. The Poisson-Fermi differential equation-that includes both steric and correlation effects-is then used to obtain the spatial profiles of energy, concentration, and dielectric coefficient from the solvent region to the filter. The energy profiles of ions are shown to depend sensitively on the steric energy that is not taken into account in the classical rate theory. We improve the rate theory by introducing a steric energy that lumps the effects of excluded volumes of all ions and water molecules and empty spaces between particles created by Lennard-Jones type and electrostatic forces. We show that the energy landscape varies significantly with bath concentrations. The energy landscape is not constant.

Ion permeation in the NanC porin from Escherichia coli: free energy calculations along pathways identified by coarse-grain simulations.

Using the X-ray structure of a recently discovered bacterial protein, the N-acetylneuraminic acid-inducible channel (NanC), we investigate computationally K(+) and Cl(-) ions' permeation. We identify ion permeation pathways that are likely to be populated using coarse-grain Monte Carlo simulations. Next, we use these pathways as reaction coordinates for umbrella sampling-based free energy simulations. We find distinct tubelike pathways connecting specific binding sites for K(+) and, more pronounced, for Cl(-) ions. Both ions permeate the porin preserving almost all of their first hydration shell. The calculated free energy barriers are G(#) ≈ 4 kJ/mol and G(#) ≈ 8 kJ/mol for Cl(-) and K(+), respectively. Within the approximations associated with these values, discussed in detail in this work, we suggest that the porin is slightly selective for Cl(-) versus K(+). Our suggestion is consistent with the experimentally observed weak Cl(-) over K(+) selectivity. A rationale for the latter is suggested by a comparison with previous calculations on strongly anion selective porins.

Correlated ions in a calcium channel model: a Poisson-Fermi theory.

We derive a continuum model, called the Poisson-Fermi equation, of biological calcium channels (of cardiac muscle, for example) designed to deal with crowded systems in which ionic species and side chains nearly fill space. The model is evaluated in three dimensions. It includes steric and correlation effects and is derived from classical hard-sphere lattice models of configurational entropy of finite size ions and solvent molecules. The maximum allowable close packing (saturation) condition is satisfied by all ionic species with different sizes and valences in a channel system, as shown theoretically and numerically. Unphysical overcrowding ("divergence") predicted by the Gouy-Chapman diffuse model (produced by a Boltzmann distribution of point charges with large potentials) does not occur with the Fermi-like distribution. Using probability theory, we also provide an analytical description of the implicit dielectric model of ionic solutions that gives rise to a global and a local formula for the chemical potential. In this primitive model, ions are treated as hard spheres and solvent molecules are described as a dielectric medium. The Poisson-Fermi equation is a local formula dealing with different correlations at different places. The correlation effects are apparent in our numerical results. Our results show variations of dielectric permittivity from bath to channel pore described by a new dielectric function derived as an output from the Poisson-Fermi analysis. The results are consistent with existing theoretical and experimental results. The binding curve of Poisson-Fermi is shown to match Monte Carlo data and illustrates the anomalous mole fraction effect of calcium channels, an effective blockage of permeation of sodium ions by a tiny concentration (or number) of calcium ions.

A parallel finite element simulator for ion transport through three-dimensional ion channel systems.

A parallel finite element simulator, ichannel, is developed for ion transport through three-dimensional ion channel systems that consist of protein and membrane. The coordinates of heavy atoms of the protein are taken from the Protein Data Bank and the membrane is represented as a slab. The simulator contains two components: a parallel adaptive finite element solver for a set of Poisson-Nernst-Planck (PNP) equations that describe the electrodiffusion process of ion transport, and a mesh generation tool chain for ion channel systems, which is an essential component for the finite element computations. The finite element method has advantages in modeling irregular geometries and complex boundary conditions. We have built a tool chain to get the surface and volume mesh for ion channel systems, which consists of a set of mesh generation tools. The adaptive finite element solver in our simulator is implemented using the parallel adaptive finite element package Parallel Hierarchical Grid (PHG) developed by one of the authors, which provides the capability of doing large scale parallel computations with high parallel efficiency and the flexibility of choosing high order elements to achieve high order accuracy. The simulator is applied to a real transmembrane protein, the gramicidin A (gA) channel protein, to calculate the electrostatic potential, ion concentrations and I - V curve, with which both primitive and transformed PNP equations are studied and their numerical performances are compared. To further validate the method, we also apply the simulator to two other ion channel systems, the voltage dependent anion channel (VDAC) and α-Hemolysin (α-HL). The simulation results agree well with Brownian dynamics (BD) simulation results and experimental results. Moreover, because ionic finite size effects can be included in PNP model now, we also perform simulations using a size-modified PNP (SMPNP) model on VDAC and α-HL. It is shown that the size effects in SMPNP can effectively lead to reduced current in the channel, and the results are closer to BD simulation results.

Ionic interactions in biological and physical systems: a variational treatment.

Chemistry is about chemical reactions. Chemistry is about electrons changing their configurations as atoms and molecules react. Chemistry has for more than a century studied reactions as if they occurred in ideal conditions of infinitely dilute solutions. But most reactions occur in salt solutions that are not ideal. In those solutions everything (charged) interacts with everything else (charged) through the electric field, which is short and long range extending to the boundaries of the system. Mathematics has recently been developed to deal with interacting systems of this sort. The variational theory of complex fluids has spawned the theory of liquid crystals (or vice versa). In my view, ionic solutions should be viewed as complex fluids, particularly in the biological and engineering context. In both biology and electrochemistry ionic solutions are mixtures highly concentrated (to approximately 10 M) where they are most important, near electrodes, nucleic ids, proteins, active sites of enzymes, and ionic channels. Ca2+ is always involved in biological solutions because the concentration (really free energy per mole) of Ca2+ in a particular location is the signal that controls many biological functions. Such interacting systems are not simple fluids, and it is no wonder that analysis of interactions, such as the Hofmeister series, rooted in that tradition has not succeeded as one would hope. Here, we present a variational treatment of ard spheres in a frictional dielectric with the hope that such a treatment of an lectrolyte as a complex fluid will be productive. The theory automatically extends to spatially nonuniform boundary conditions and the nonequilibrium systems and flows they produce. The theory is unavoidably self-consistent since differential equations are derived (not assumed) from models of (Helmholtz free) nergy and dissipation of the electrolyte. The origin of the Hofmeister series is (in my view) an inverse problem that becomes well posed when enough data from disjoint experimental traditions are interpreted with a self-consistent theory.

Interacting ions in biophysics: real is not ideal.

Ions in water are important throughout biology, from molecules to organs. Classically, ions in water were treated as ideal noninteracting particles in a perfect gas. Excess free energy of each ion was zero. Mathematics was not available to deal consistently with flows, or interactions with other ions or boundaries. Nonclassical approaches are needed because ions in biological conditions flow and interact. The concentration gradient of one ion can drive the flow of another, even in a bulk solution. A variational multiscale approach is needed to deal with interactions and flow. The recently developed energetic variational approach to dissipative systems allows mathematically consistent treatment of the bio-ions Na(+), K(+), Ca(2+), and Cl(-) as they interact and flow. Interactions produce large excess free energy that dominate the properties of the high concentration of ions in and near protein active sites, ion channels, and nucleic acids: the number density of ions is often >10 M. Ions in such crowded quarters interact strongly with each other as well as with the surrounding protein. Nonideal behavior found in many experiments has classically been ascribed to allosteric interactions mediated by the protein and its conformation changes. The ion-ion interactions present in crowded solutions-independent of conformation changes of the protein-are likely to change the interpretation of many allosteric phenomena. Computation of all atoms is a popular alternative to the multiscale approach. Such computations involve formidable challenges. Biological systems exist on very different scales from atomic motion. Biological systems exist in ionic mixtures (like extracellular and intracellular solutions), and usually involve flow and trace concentrations of messenger ions (e.g., 10(-7) M Ca(2+)). Energetic variational methods can deal with these characteristic properties of biological systems as we await the maturation and calibration of all-atom simulations of ionic mixtures and divalents.

Ionic interactions are everywhere.

Ionic solutions are dominated by interactions because they must be electrically neutral, but classical theory assumes no interactions. Biological solutions are rather like seawater, concentrated enough so that the diameter of ions also produces important interactions. In my view, the theory of complex fluids is needed to deal with the interacting reality of biological solutions.

PNP equations with steric effects: a model of ion flow through channels.

The flow of current through an ionic channel is studied using the energetic variational approach of Liu applied to the primitive (implicit solvent) model of ionic solutions. This approach allows the derivation of self-consistent (Euler-Lagrange) equations to describe the flow of spheres through channels. The partial differential equations derived involve the global interactions of the spheres and are replaced here with a local approximation that we call steric PNP (Poisson-Nernst-Planck) (Lin, T. C.; Eisenberg, B. To be submitted for publication, 2012). Kong combining rules are used and a range of values of steric interaction parameters are studied. These parameters change the energetics of steric interaction but have no effect on diffusion coefficients in models and simulations. Calculations are made for the calcium (EEEE, EEEA) and sodium channels (DEKA) previously studied in Monte Carlo simulations with comparable results. The biological function is quite sensitive to the steric interaction parameters, and we speculate that a wide range of the function of channels and transporters, even enzymes, might depend on such terms. We point out that classical theories of channels, transporters, and enzymes depend on ideal representations of ionic solutions in which nothing interacts with nothing, even in the enormous concentrations found near and in these proteins or near electrodes in electrochemical cells for that matter. We suggest that a theory designed to handle interactions might be more appropriate. We show that one such theory is feasible and computable: steric PNP allows a direct comparison with experiments measuring flows as well as equilibrium properties. Steric PNP combines atomic and macroscales in a computable formulation that allows the calculation of the macroscopic effects of changes in atomic scale structures (size ~/= 10(-10) meters) studied very extensively in channology and molecular biology.

Ionizable side chains at catalytic active sites of enzymes.

Catalytic active sites of enzymes of known structure can be well defined by a modern program of computational geometry. The CASTp program was used to define and measure the volume of the catalytic active sites of 573 enzymes in the Catalytic Site Atlas database. The active sites are identified as catalytic because the amino acids they contain are known to participate in the chemical reaction catalyzed by the enzyme. Acid and base side chains are reliable markers of catalytic active sites. The catalytic active sites have 4 acid and 5 base side chains, in an average volume of 1,072 Å(3). The number density of acid side chains is 8.3 M (in chemical units); the number density of basic side chains is 10.6 M. The catalytic active site of these enzymes is an unusual electrostatic and steric environment in which side chains and reactants are crowded together in a mixture more like an ionic liquid than an ideal infinitely dilute solution. The electrostatics and crowding of reactants and side chains seems likely to be important for catalytic function. In three types of analogous ion channels, simulation of crowded charges accounts for the main properties of selectivity measured in a wide range of solutions and concentrations. It seems wise to use mathematics designed to study interacting complex fluids when making models of the catalytic active sites of enzymes.